#include "ellipse.h"

#define NUM_SEGMENTS 32
#define W 32  // 使用 32 位处

#define CRYPTO_START_SEC_CODE
#include "Crypto_MemMap.h"

/**
* Parameter Constant
*/

// 椭圆曲线参数
const sm2_uint256_t Prime = {
    0xFFFFFFFE, 0xFFFFFFFF,
    0xFFFFFFFF, 0xFFFFFFFF,
    0xFFFFFFFF, 0x00000000,
    0xFFFFFFFF, 0xFFFFFFFF
};
const sm2_uint256_64t Prime_64 = {
    0xFFFFFFFEFFFFFFFF,
    0xFFFFFFFFFFFFFFFF,
    0xFFFFFFFF00000000,
    0xFFFFFFFFFFFFFFFF
};
const uint8_t param_A[32] = {
                0xFF,0xFF,0xFF,0xFE,0xFF,0xFF,0xFF,0xFF,0xFF,0xFF,0xFF,0xFF,0xFF,0xFF,0xFF,0xFF,
                0xFF,0xFF,0xFF,0xFF,0x00,0x00,0x00,0x00,0xFF,0xFF,0xFF,0xFF,0xFF,0xFF,0xFF,0xFC, };
const uint32_t _32b_param_A[8] = {
                0xFFFFFFFE,0xFFFFFFFF,0xFFFFFFFF,0xFFFFFFFF,
                0xFFFFFFFF,0x00000000,0xFFFFFFFF,0xFFFFFFFC, };
const sm2_uint256_64t param_A_64 = {
                0xFFFFFFFEFFFFFFFF,0xFFFFFFFFFFFFFFFF,
                0xFFFFFFFF00000000,0xFFFFFFFFFFFFFFFC, };

const uint32_t mg_A[8] = {
    0xfffffffb, 0xffffffff, 0xffffffff, 0xffffffff,
    0xfffffffc, 0x00000003, 0xffffffff, 0xfffffffc };
const sm2_uint256_64t mg_A_64 = {
    0xfffffffbffffffff, 0xffffffffffffffff,
    0xfffffffc00000003, 0xfffffffffffffffc };

const uint8_t param_B[32] = {
                0x28,0xE9,0xFA,0x9E,0x9D,0x9F,0x5E,0x34,0x4D,0x5A,0x9E,0x4B,0xCF,0x65,0x09,0xA7,
                0xF3,0x97,0x89,0xF5,0x15,0xAB,0x8F,0x92,0xDD,0xBC,0xBD,0x41,0x4D,0x94,0x0E,0x93, };
const sm2_uint256_64t mg_B_64 = {
    0x240fe188ba20e2c8, 0x527981505ea51c3c, 0x71cf379ae9b537ab, 0x90d230632bc0dd42,
};

const sm2_point_64t param_G_64 = {
    .x = {0x32C4AE2C1F198119, 0x5F9904466A39C994,
          0x8FE30BBFF2660BE1, 0x715A4589334C74C7,},
    .y = {0xBC3736A2F4F6779C, 0x59BDCEE36B692153,
          0xD0A9877CC62A4740, 0x02DF32E52139F0A0,},
    .z = {0,0,0,1}
};

const sm2_point_64t mg_G_64 = {
    .x = {0x91167a5ee1c13b05, 0xd6a1ed99ac24c3c3,
          0x3e7981eddca6c050, 0x61328990f418029e},
    .y = {0x63cd65d481d735bd, 0x8d4cfb066e2a48f8,
          0xc1f5e5788d3295fa, 0xc1354e593c2d0ddd},
    .z = {
         0x0000000100000000,0000000000000000,0x00000000ffffffff,0000000000000001,
    }
};

const sm2_uint256_t param_N = {
            0xFFFFFFFE, 0xFFFFFFFF, 0xFFFFFFFF, 0xFFFFFFFF,
            0x7203DF6B, 0x21C6052B, 0x53BBF409, 0x39D54123
};
const sm2_uint256_64t param_N_64 = {
            0xFFFFFFFEFFFFFFFF, 0xFFFFFFFFFFFFFFFF,
            0x7203DF6B21C6052B, 0x53BBF40939D54123
};


const sm2_uint256_64t P_NEG_INV_64 = {
    0xfffffffc00000001,
    0xfffffffe00000000,
    0xffffffff00000001,
    0x0000000000000001
};

const uint64_t N_NEG_INV_64[4] = {
    0x6f39132f82e4c7bc,
    0x2b0068d3b08941d4,
    0xdf1e8d34fc8319a5,
    0x327f9e8872350975

};

const uint64_t N_2_512_mod_n_64[4] = {
    0x1eb5e412a22b3d3b,  0x620fc84c3affe0d4, 0x3464504ade6fa2fa, 0x901192af7c114f20
};

const uint64_t N_2_512_mod_p_64[4] = {
    0x0000000400000002,  0x0000000100000001, 0x00000002ffffffff, 0x0000000200000003
};

const sm2_point_64t ZERO_POINT_64 = {
    .x = {0,0,0,0},
    .y = {0,0,0,0},
    .z = {0,0,0,0},
};

// 引入256位数做计算  数据存储为大端存储  高位---------低位

const sm2_uint256_64t N_256_1_64 = { 0,0,0,1 };
const sm2_uint256_64t N_256_ZERO_64 = { 0,0,0,0 };

sm2_uint256_64t num_3_64 = { 0x0000000300000000,0x0000000000000000,0x00000002fffffffd,0x0000000000000003 };
sm2_uint256_64t num_2_64 = { 0x0000000200000000,0x0000000000000000,0x00000001fffffffe,0x0000000000000002 };
sm2_uint256_64t N_ONE_64 = { 0x0000000100000000,0x0000000000000000,0x00000000ffffffff,0x0000000000000001 };
sm2_uint256_64t N_TWO_INV_64 = { 0x8000000000000000, 0x0000000000000000, 0x0000000000000000, 0x0000000000000000 };
sm2_uint256_64t N_ONE_MOD_n_64 = { 0x0000000100000000, 0x0000000000000000, 0x8dfc2094de39fad4, 0xac440bf6c62abedd, };
// 多倍点预计算
// 
sm2_point_64t precomputed_table[17] = {
{
    .x = {0x0000000000000000, 0x0000000000000000, 0x0000000000000000, 0x0000000000000000,
},

. y = {0x0000000000000000, 0x0000000000000000, 0x0000000000000000, 0x0000000000000000,
},

. z = {0x0000000000000000, 0x0000000000000000, 0x0000000000000000, 0x0000000000000000,
},
},
{
    .x = {0x91167a5ee1c13b05, 0xd6a1ed99ac24c3c3, 0x3e7981eddca6c050, 0x61328990f418029e,
},

. y = {0x63cd65d481d735bd, 0x8d4cfb066e2a48f8, 0xc1f5e5788d3295fa, 0xc1354e593c2d0ddd,
},

. z = {0x0000000100000000, 0x0000000000000000, 0x00000000ffffffff, 0x0000000000000001,
},
},
{
    .x = {0x398874c476a3b1f7, 0x7aef3e8626014409, 0x03243d78d5b614a6, 0x2eda8381e63c48d6,
},

. y = {0x1fbbdfdddaf4fd47, 0x5a86a7ae64921d48, 0x29f04a88f6cf4dc1, 0x28385681c1a73e40,
},

. z = {0xc79acba903ae6b7b, 0x1a99f60cdc5491f1, 0x83ebcaf11a652bf5, 0x826a9cb2785a1bba,
},
},
{
    .x = {0x3015ad2077c8c028, 0x5dd76fbcbdf037e4, 0x3c51eb40c7c50abc, 0x05294c0f53498d0e,
},

. y = {0x167df558de1a1105, 0x86371f717ce702b4, 0x9efa98a8c6242c3c, 0x0d40bb897bf1733c,
},

. z = {0xe3d2125af942174e, 0x02c20fad90368de2, 0x5832932b7eef2b43, 0x51711cc32c29b331,
},
},
{
    .x = {0x87b2ca9ded2487c6, 0xefdbc69303258763, 0xa0b5520fc63cf401, 0x54f6c059b945acf2,
},

. y = {0xdc86353bc72db45e, 0xbb5b2d03cec4614b, 0x164688f19f12dd85, 0x7fd007e181457b59,
},

. z = {0x050653f8579d1d2d, 0x930d7346e31bad56, 0xb5a4654d6a9f2c50, 0x22434941744ced3a,
},
},
{
    .x = {0x6cd1f0d1f755d65e, 0x2b5b7201cc2acb83, 0x3c5884751bb71452, 0xf6ddc25c1dbebcb3,
},

. y = {0x1424f18af712f779, 0x4329f78af9786e24, 0x15cb971ea9880e78, 0xb21ce1c564ceb2d4,
},

. z = {0x76e2d7521adac9a9, 0xcf768feaf4f3c22b, 0x042bb8f4e5f88337, 0xb3fa76369f25508f,
},
},
{
    .x = {0xb72690c01db7d188, 0xa107b71dfa3a9318, 0x297e8f1473a4996b, 0x22f45b47df293cd2,
},

. y = {0x585770c411ee9959, 0xacc70f7886206512, 0x16d8522374354144, 0x978fc7d13515abb5,
},

. z = {0x65a472cea9f5c346, 0x599e767846d6d6a9, 0xcbde61278cf02456, 0x9b05a718655f9aa4,
},
},
{
    .x = {0xc438f40566154ee6, 0x092e1e43f2aadc29, 0xd1f4bfff9d80615f, 0xfbebd9b221a1dc21,
},

. y = {0x48507ea2be69f91f, 0x8dea43f44df05372, 0xea675e9186a70704, 0xa5c8bf390113f620,
},

. z = {0xa367a69298acadf2, 0xd911ad3528a331ac, 0xcec0b9f8ca50ebba, 0x164c8e57cf549bec,
},
},
{
    .x = {0xf5c816a1a84e6d13, 0xdf9ad5506b788970, 0xc02ae0c983a0e64d, 0x0590b9d3af0989e6,
},

. y = {0x2cf8d1491bafee92, 0x917a4d4651005d24, 0x14668d19bdb94608, 0x88841f170d7d7614,
},

. z = {0xa3325160280d774e, 0x4ec7df5751da58da, 0xbecd5bfdaa66ad8b, 0x249f6a348953a88a,
},
},
{
    .x = {0x224fabd16ea55eec, 0x6e44e0a36e36b39a, 0xeb03fc82868ddd04, 0x006f81f8bedf1939,
},

. y = {0xadafb0fa363d2328, 0x63d5b0f20f45dbe4, 0x4010c2a2c0d60fe0, 0x191e56b17030b4be,
},

. z = {0xffe14fe053771949, 0xea502c884e3b254c, 0x493c4aa64e28ce5d, 0xda66802983580a49,
},
},
{
    .x = {0xdb312d34f7f0a8b8, 0xebcc7980ee139f89, 0x537940e8f48d8054, 0xa9831c06366c16ec,
},

. y = {0x622d9c599d9d31c5, 0x8d6c1fc624abf9ed, 0x3ad88554ebaca6c4, 0x7802ceea0312ccbf,
},

. z = {0xf31fa36b57da7db2, 0x1a268162f64e9178, 0xb2c395c8eeb668ca, 0x70c0d7a50d54e8c6,
},
},
{
    .x = {0xcf9c805da845b996, 0x07a2a57bcb0f24f6, 0xa5352f81d6cf3f46, 0x226bd219cd6f9614,
},

. y = {0x13eceefa0c1073e4, 0x632eee5964b0aa84, 0x47bddf0a11e4289c, 0x2a6a4d1fe37b4c58,
},

. z = {0x78676a3c0a581003, 0x78c4298242c4df85, 0x8249b9e07e1fd290, 0xf66d84677424943b,
},
},
{
    .x = {0xb14233809afee301, 0x338b60f2cf1c1dc0, 0x05f6f0e36320c4d0, 0x32891b12c905a575,
},

. y = {0xddf9f0c402d7f8d3, 0x61ca440926ce0b7c, 0x155f095057df498c, 0x9ccb3884bd05288c,
},

. z = {0x7d8b33476859c5f7, 0xbafb675055b64f07, 0x2cf17c4518a4fec0, 0x0f01014f83fa1805,
},
},
{
    .x = {0xee51577374602439, 0xc88124ae212f8c56, 0x6e194f729e1905ae, 0xd8cdeb32279f5515,
},

. y = {0x7e47adce2abfd7fa, 0x1e1eb4e7f65ab8c4, 0x13f16ecd9030b9d6, 0x8fa84c7282cc4b08,
},

. z = {0xf9e66437dbccf742, 0xdf2ba81113c5dc37, 0xa976d0ec9ac4f0fd, 0x44fd7942fffd1a28,
},
},
{
    .x = {0x002cd2339238095c, 0x851747e920d28121, 0xa348e6cf61d63281, 0x375e242a3185b621,
},

. y = {0xbc32ac1ab4b09567, 0x3f1f101560a893ca, 0xbed4fe9257db3c46, 0xbdad5a036794241e,
},

. z = {0x1765d20a1eeb66d1, 0xb8a0cae15406d7ae, 0xf80d47b33ebf03c1, 0xecf6743f506d5e2e,
},
},
{
    .x = {0x6fb69aaf5095848c, 0xc5d70a301b825aca, 0xd6f0b0d2dd7d9c2e, 0x0810b179ec4f486b,
},

. y = {0x25e2b07654be2acc, 0x49d94ab82fd141b1, 0xc251e10839277156, 0x53c1b3b37837f8fc,
},

. z = {0x0c8d9435af5f5a29, 0xc3185b0f82c8feba, 0xbfd04afade432751, 0xcc0a55e63d8c2e89,
},
},
{
    .x = {0x54ea747ff266e11d, 0x75c5d1db1c5abe01, 0x953dc915f0df5af3, 0xe899b82f1a71c36f,
},

. y = {0x211f9870979316c0, 0x044d1bea3613f6fa, 0x0f5c0a0b1629eacd, 0x7ce89b13b22caba5,
},

. z = {0x1acc5cec17e31a58, 0x1de8319dbd637ee5, 0xe466ae787547e09f, 0xa5a2c36b856379fa,
},
},
};
void rightShift256(uint64_t* bigEndian256, int shiftAmount) {
    uint64_t carry = 0;
    for (int i = 0; i <4; i++) {
        uint64_t newCarry = (bigEndian256[i] << (64 - shiftAmount)) & 0xFFFFFFFFFFFFFFFF;
        bigEndian256[i] = (bigEndian256[i] >> shiftAmount) | carry;
        carry = newCarry;
    }
}

 //基于 wNAF（Window Non - Adjacent Form）优化的标量乘法实现
sm2_point_64t  Crypto_ECC_wNAF_scalar_multiply(sm2_point_64t P, sm2_uint256_64t k, int w) {
    int N = 1 << w;  // 2^w，窗口大�?
    int naf[256] = { 0 }; // 存储 wNAF 的表�?
    //sm2_point_t* precomputed = (sm2_point_t*)malloc(N * sizeof(sm2_point_t));  // 存储预计算的倍点
    //sm2_point_64t precomputed_table[256];

    //// Step 1: 计算 wNAF
    uint64_t index = 0;
    uint64_t temp_k, carry=0;
    sm2_uint256_64t t_k;
    (void)ILib_memcpy(t_k, k, sizeof(sm2_uint256_64t));
    
    temp_k = k[0];
    carry = 0;
    while (t_k[3]>0) {
        if (t_k[3] & 1) {
            int e = t_k[3] & ((1 << w)-1);  // 计算k�?2^w
            if (e >= (1 << (w - 1))) {  // 如果e >= 2^(w-1)，则将其转换为负�?
                e -= (1 << w);  // 减去2^w得到负�?
            }
            naf[index] = e; // 非邻接的形式�?1表示该位�?1
            //temp_k = temp_k - naf[index];
            t_k[3] -= e;
        }
        else {
            naf[index] = 0;
        }
        index++;
        rightShift256(t_k, 1);
        //temp_k >>= 1;
            
    }

    //printf("mod %u: [", index);
    //for (int j = index - 1; j >= 0; j--) {
    //    printf("%d", naf[j]);
    //    if (j > 0) {
    //        printf(", ");
    //    }
    //}
    //printf("]\n");

    // Step 3: 执行 wNAF 标量乘法
    //sm2_point_64t R = (sm2_point_64t){ {0} };  // 初始化为零点
    sm2_point_64t Q;
    (void)ILib_memset(&Q, 0, sizeof(sm2_point_64t));
    sm2_point_64t temp;
    for (int i = index-1; i >=0; i--) {
        // 执行�?2操作：两倍点
        Crypto_ECC_jacc_point_double_64(&Q, Prime_64, &Q);
        //Q = precomputed_sqr_table[cnt];
        //cnt++;
        int value = naf[i];
        if (value != 0) {
            if (value > 0) {
                Crypto_ECC_jacc_point_add_64(&Q, &precomputed_table[abs(value)], Prime_64, &Q);  // 加上对应的预计算倍点
            }
            if (value < 0) {
                point_inverse_64(&precomputed_table[abs(value)], Prime_64, &temp);
                //print_point(&temp);
                Crypto_ECC_jacc_point_add_64(&Q, &temp, Prime_64, &Q);  // 加上对应的预计算倍点
            }
        }
    }

    return Q;
}

/**
 * Convert the uint8_t array back to the uint32_t array
 *
 * @param src uint32_t array to be converted
 * @param dst uint8_t array of the conversion
 */
void convert_to_uint32(uint8_t src[32], uint32_t dst[8]) {
    for (size_t i = 0; i < 8; i++) {
        uint32_t value = 0;
        for (size_t j = 0; j < 4; j++) {
            value = (value << 8) | src[i * 4 + j];
        }
        dst[i] = value;
    }
}
void convert_uint8_to_uint64(uint8_t src[32], uint64_t dst[4]) {
    for (size_t i = 0; i < 4; i++) {
        uint64_t value = 0;
        for (size_t j = 0; j < 8; j++) {
            value = (value << 8) | src[i * 8 + j];
        }
        dst[i] = value;
    }
}
void convert_uint64_to_uint8(uint64_t src[4], uint8_t dst[32]) {
    for (size_t i = 0; i < 4; i++) {
        dst[i * 8 + 0] = (src[i] >> 56) & 0xFF;
        dst[i * 8 + 1] = (src[i] >> 48) & 0xFF;
        dst[i * 8 + 2] = (src[i] >> 40) & 0xFF;
        dst[i * 8 + 3] = (src[i] >> 32) & 0xFF;
        dst[i * 8 + 4] = (src[i] >> 24) & 0xFF;
        dst[i * 8 + 5] = (src[i] >> 16) & 0xFF;
        dst[i * 8 + 6] = (src[i] >> 8) & 0xFF;
        dst[i * 8 + 7] = (src[i]) & 0xFF;
    }
}

void convert_uint32_to_uint64(uint32_t src[8], uint64_t dst[4]) {
    for (size_t i = 0; i < 4; i++) {
        // 利用类型转换，直接将�?32位和�?32位拼接成一�?64位数�?
        dst[i] = ((uint64_t)src[2 * i]) << 32 | (((uint64_t)src[2 * i + 1]));
    }
}



/**
 * Convert the uint32_t array back to the uint8_t array
 *
 * @param src uint8_t array to be converted
 * @param dst uint32_t array of the conversion
 */
void convert_to_uint8(uint32_t src[8], uint8_t dst[32]) {

    for (size_t i = 0; i < 8; i++) {
        for (size_t j = 0; j < 4; j++) {
            dst[i * 4 + j] = (src[i] >> (24 - j * 8)) & 0xFF;
        }
    }
}


 /**
   * Convert number from the natural domain to number in the Montgomery domain(P)
   *
   *
   * @param src uint32_t array in natural domain to be converted
   * @param p prime which needed from Montgomery domain
   * @param dst uint32_t array in Montgomery domain of the conversion
   */

void num2mgmy_modp_64(uint64_t src[4], sm2_uint256_64t p, uint64_t dst[4])
{
    //uint32_t N_2_512[8];
    montgomery_multiplication_modp_64_32(src, N_2_512_mod_p_64, p, dst);


}
void num2mgmy_modn_64(uint64_t src[4], sm2_uint256_64t p, uint64_t dst[4])
{
     montgomery_multiplication_modn_64_32(src, N_2_512_mod_n_64, p, dst);
}
/**
 * Convert number from the Montgomery domain(N) to number in the natural domain
 *
 * @param src uint32_t array in Montgomery domain(N) to be converted
 * @param p prime which needed from Montgomery domain
 * @param dst uint32_t array in natural domain of the conversion
 */

void mgmy2num_modp_64(uint64_t src[4], sm2_uint256_64t p, uint64_t dst[4])
{
    montgomery_multiplication_modp_64_32(N_256_1_64, src, p, dst);

    //mod_mul(src, N_256_1, p, dst);
}
void mgmy2num_modn_64(uint64_t src[4], sm2_uint256_64t p, uint64_t dst[4])
{
    montgomery_multiplication_modn_64_32(N_256_1_64, src, p, dst);
    //mod_mul(src, N_256_1, p, dst);
}
/**
 * Convert point from the natural domain to number in the Montgomery domain
 *
 *
 * @param src point in natural domain to be converted
 * @param p prime which needed from Montgomery domain
 * @param dst point in Montgomery domain of the conversion
 */


void point_modp2mgmy_64(sm2_point_64t* src, sm2_uint256_64t p, sm2_point_64t* dst) {
    if (src == NULL || dst == NULL) {
        return;
    }
    num2mgmy_modp_64(src->x, p, dst->x);
    num2mgmy_modp_64(src->y, p, dst->y);
    (void)ILib_memcpy(dst->z, N_ONE_64, sizeof(sm2_uint256_64t));
}

/**
 * Convert point from the Montgomery domain to number in the natural domain
 *
 * @param src point in Montgomery domain to be converted
 * @param p prime which needed from Montgomery domain
 * @param dst point in natural domain of the conversion
 */

void mgmy_modp2point_64(sm2_point_64t* src, sm2_uint256_64t p, sm2_point_64t* dst) {
    if (src == NULL || dst == NULL) {
        return;
    }
    Crypto_ECC_jacobi_to_affine_64(src->x, src->y, src->z, p, src);
    mgmy2num_modp_64(src->x, p, dst->x);
    mgmy2num_modp_64(src->y, p, dst->y);
}

/**
 * Take the lower 256 bits of 512 digits
 *
 * @param src 512 digits
 * @param dst Low 256 bits of src
 */

void sm2_getlow256_64(sm2_uint512_64t src, sm2_uint256_64t dst) {
    if (src == NULL || dst == NULL) {
        return;
    }
    for (int i = 4; i < 8; i++)
    {
        dst[i - 4] = src[i];
    }
}

/**
 * Take the higher 256 bits of 512 digits
 *
 * @param src 512 digits
 * @param dst high 256 bits of src
 */

void sm2_gethigh256_64(sm2_uint512_64t src, sm2_uint256_64t dst) {
    if (src == NULL || dst == NULL) {
        return;
    }
    for (int i = 0; i < 4; i++)
    {
        dst[i] = src[i];
    }
}
/**
 * Compare the size of 256 bits, compare every 32 bits. 
 *
 * @Return  a>b return 1；a<b return -1；a=b return 0
 */

int compare_value(sm2_uint256_t a, sm2_uint256_t b) {
    for (int i = 0; i < 8; i++)
    {
        if (a[i] > b[i])
            return 1;
        else if (a[i] < b[i])
            return -1;
    }
    return 0;
}
int compare_value_64(sm2_uint256_64t a, sm2_uint256_64t b) {
    for (int i = 0; i < 4; i++)
    {
        if (a[i] > b[i])
            return 1;
        else if (a[i] < b[i])
            return -1;
    }
    return 0;
}


 /**
  * Invert each element of P
  */
void neg_p_64(sm2_uint256_64t result, const sm2_uint256_64t p) {
    // �? P 每个元素进行取反
    for (int i = 0; i < 4; i++) {
        result[i] = ~p[i];
    }

    // �? 1 处理，从最低有效位开�?
    uint64_t carry = 1;
    for (int i = 3; i >= 0 && carry != 0; i--) {
        uint64_t sum = result[i] + carry;
        result[i] = sum;
        carry = (sum < carry);  // 如果发生溢出，则继续进位
    }
}

/**
 * Addition and modulo operation of two numbers
 *
 * @param a addend number 1
 * @param b addend number 2
 * @param p prime
 * @param result the result of adding two numbers and taking the modulus
 */
uint64_t add_64(sm2_uint256_64t a, sm2_uint256_64t b, sm2_uint256_64t result)
{
    uint64_t carry = 0;  // 初始化进�?
    for (int i = 3; i >= 0; i--) {
        uint64_t sum = a[i] + b[i] + carry;  // 当前位的加法
        carry = (sum < a[i]) || (sum < b[i]) || (sum == a[i] && carry) || (sum == b[i] && carry);  // 检测是否进�?
        result[i] = sum;  // 存储结果
    }
    return carry;
}
void mod_add_64(sm2_uint256_64t a, sm2_uint256_64t b, sm2_uint256_64t p, sm2_uint256_64t result) {

    uint64_t carry = add_64(a, b, result);

    // 进位或结果超出模数，直接减去模数
    if (carry != 0 || compare_value_64(result, p) >= 0) {
        sub_64(result, p, result);
       }
}

/**
 * Ordinary addition
 */

uint64_t add_128(sm2_uint128_t a, sm2_uint128_t b, sm2_uint128_t result)
{
    uint64_t carry = 0;  // 初始化进�?
    for (int i = 1; i >= 0; i--) {
        uint64_t sum = a[i] + b[i] + carry;  // 当前位的加法
        //carry = (sum < a[i]) || (sum < b[i]) || (sum==a[i]&&carry)|| (sum == b[i] && carry); // 检测是否进�?
        carry = (sum < a[i]) || (sum < b[i]);
        result[i] = sum;  // 存储结果
    }
    return carry;
}

/**
 * 512 digit addition operation
 */
uint32_t add_uint512_t(sm2_uint512_t a, sm2_uint512_t b, sm2_uint512_t result)
{
    uint32_t carry = 0;  // 初始化进�?
    for (int i = 15; i >= 0; i--) {
        uint32_t sum = a[i] + b[i] + carry;  // 当前位的加法
        carry = (sum < a[i]) || (sum < b[i] && carry);  // 检测是否进�?
        result[i] = sum;  // 存储结果
    }
    return carry;
}

uint64_t add_uint512_64t(sm2_uint512_64t a, sm2_uint512_64t b, sm2_uint512_64t result)
{
    uint64_t carry = 0;  // 初始化进�?
    for (int i = 7; i >= 0; i--) {
        uint64_t sum = a[i] + b[i] + carry;  // 当前位的加法
        carry = (sum < a[i]) || (sum < b[i] && carry);  // 检测是否进�?
        result[i] = sum;  // 存储结果
    }
    return carry;
}
/**
 * Subtraction and modulo operation of two numbers
 *
 * @param a minuend number 
 * @param b subtrahend number 
 * @param p prime
 * @param result the result of subtracting two numbers and taking the modulus
 */

void mod_sub_64(sm2_uint256_64t a, sm2_uint256_64t b, sm2_uint256_64t p, sm2_uint256_64t result)
{
    uint64_t borrow = sub_64(a, b, result);
    // 有借位
    if (borrow)
    {
        sm2_uint256_64t np;
        neg_p_64(np, p);
        sub_64(result, np, result);
    }
}

/**
 * Ordinary subtraction
 */
uint32_t sub(sm2_uint256_t a, sm2_uint256_t b, sm2_uint256_t result)
{
    sm2_uint256_t res;
    uint32_t borrow = 0;
    uint32_t tp;
    //memset(res, 0, sizeof(sm2_uint256_t));
    (void)ILib_memset(res, 0, sizeof(sm2_uint256_t));
    tp = a[7] - b[7];
    res[7] = tp;
    borrow = (tp > a[7]);
    for (int i = 6; i >= 0; i--) {
        tp = a[i] - borrow;
        borrow = (tp > a[i]);
        res[i] = tp - b[i];
        borrow += (res[i] > tp); // 检查是否需要借位
    }
    //memcpy(result, res, sizeof(sm2_uint256_t));
    ILib_memcpy(result, res, sizeof(sm2_uint256_t));
    return borrow;
}
uint64_t sub_64(sm2_uint256_64t a, sm2_uint256_64t b, sm2_uint256_64t result)
{
    sm2_uint256_64t res;
    uint64_t borrow = 0;
    uint64_t tp;
    //memset(res, 0, sizeof(sm2_uint256_t));
    (void)ILib_memset(res, 0, sizeof(sm2_uint256_64t));
    tp = a[3] - b[3];
    res[3] = tp;
    borrow = (tp > a[3]);
    for (int i = 2; i >= 0; i--) {
        tp = a[i] - borrow;
        borrow = (tp > a[i]);
        res[i] = tp - b[i];
        borrow += (res[i] > tp); // 检查是否需要借位
    }
    //memcpy(result, res, sizeof(sm2_uint256_t));
    (void)ILib_memcpy(result, res, sizeof(sm2_uint256_64t));
    return borrow;
}

uint32_t sub_128(sm2_uint128_t a, sm2_uint128_t b, sm2_uint128_t result)
{
    sm2_uint128_t res;
    uint64_t borrow = 0;
    uint64_t tp;
    //memset(res, 0, sizeof(sm2_uint256_t));
    (void)ILib_memset(res, 0, sizeof(sm2_uint128_t));
    tp = a[1] - b[1];
    res[1] = tp;
    borrow = (tp > a[1]);
    tp = a[0] - borrow;
    borrow = (tp > a[0]);
    res[0] = tp - b[0];
    borrow += (res[0] > tp); // 检查是否需要借位
    //memcpy(result, res, sizeof(sm2_uint256_t));
    ILib_memcpy(result, res, sizeof(sm2_uint128_t));
    return borrow;
}

/**
 * Multiplying two 256 bit numbers
 *
 * @param a 256 bit number
 * @param b 256 bit number
 * @param result multiplying two 256 bit numbers yields 512 bit numbers
 */
uint32_t mul(sm2_uint256_t a, sm2_uint256_t b, sm2_uint512_t result)
{
    uint32_t temp_result[16] = { 0 };  // 存储 512 位乘积结�?
    uint16_t* tr16 = (uint16_t*)temp_result;
    uint32_t carry = 0;  // 初始化进�?
    // 逐位计算乘积
    int i, j;
    for (i = 0; i < 16; i++) {
        carry = 0;
        for (j = 0; j < 16; j++) {
            uint32_t ta = ((a[(15 - i) / 2] >> (((16 - i) % 2) * 16)) & 0xFFFF);
            uint32_t tb = ((b[(15 - j) / 2] >> (((16 - j) % 2) * 16)) & 0xFFFF);
            uint32_t sum = tr16[i + j] + ta * tb + carry;
            tr16[i + j] = sum & 0xFFFF;
            carry = sum >> 16;
        }
        tr16[i + 16] = carry;
    }

    for (int t = 0; t < 16; t++)
    {
        result[15 - t] = temp_result[t];
    }

    return carry;
}
#  define DATA(a)        ((a)&0xffffffff)
#  define CARRY(a)        (((a)>>32))
void montgomery_multiplication_modp_64_32(uint64_t* A, uint64_t* B, uint64_t* p, uint64_t* R) {

    if (!ILib_memcmp(A, N_256_ZERO_64, sizeof(Prime_64)) || !ILib_memcmp(B, N_256_ZERO_64, sizeof(Prime_64))) {
        (void)ILib_memset(R, 0, sizeof(Prime_64));
        return;
    }
    //if (sm2_uint256t_is_zero(A) || sm2_uint256t_is_zero(B)) {
    //    sm2_uint256t_set_zero(R);
    //    return;
    //}
    int S = 8;

    uint32_t t[10] = { 0 };
    //uint32_t* t = (uint32_t*)malloc(18 * sizeof(uint32_t));
    //ILib_memset(t, 0, 18);
    uint32_t* t_base = t;
    uint32_t T[9] = { 0 };
    uint64_t C;
    //uint64_t C; // uint128
    //sm2_uint128_t C, res;
    uint64_t carry = 0;

    sm2_uint256_t n_temp, n_inv_temp, p_temp;
    sm2_uint256_64t neg_inv;
    uint64_t A_temp[8], B_temp[8];

    // 改为小端数据进行操作
    int cnt = 6;
    for (int i = 0; i < 4; i++) {
        // 小端存储，将64位数据拆分成2�?32位，并以小端格式存储
        A_temp[cnt] = (A[i] & 0xFFFFFFFF);        // �? 32 �?
        A_temp[cnt + 1] = ((A[i] >> 32) & 0xFFFFFFFF); // �? 32 �?

        B_temp[cnt] = (B[i] & 0xFFFFFFFF);        // �? 32 �?
        B_temp[cnt + 1] = ((B[i] >> 32) & 0xFFFFFFFF); // �? 32 �?

        n_temp[cnt] = (p[i] & 0xFFFFFFFF);        // �? 32 �?
        n_temp[cnt + 1] = ((p[i] >> 32) & 0xFFFFFFFF); // �? 32 �?

        n_inv_temp[cnt] = (P_NEG_INV_64[i] & 0xFFFFFFFF);    // �? 32 �?
        n_inv_temp[cnt + 1] = ((P_NEG_INV_64[i] >> 32) & 0xFFFFFFFF); // �? 32 �?

        cnt -= 2; // 由于每个64位数拆分为两�?32位数，所以减�?2
    }
    //print_uint32_array(A_temp);
    //printf("\n");
    // 外层循环
    for (size_t i = 0; i < S; i++) {
        C = 0;
        // 内层循环：计�? A[j] * B[i] 并累加到 t[j]
        for (size_t j = 0; j < S; j++) {
            C = (uint64_t)t[j] + A_temp[j] * B_temp[i] + C;
            t[j] = DATA(C);
            C >>= W;
        }

        C += (uint64_t)t[S];
        //t[S] = C & ((1 << W) - 1);
        //t[S + 1] = C >> W;
        t[S] = DATA(C);
        t[S + 1] = C >> W;
        C = 0;

        // 计算 m = t[0] * n_inv (mod 2^W)
        uint32_t m = (t[0] * n_inv_temp[0]) & 0xffffffff;

        // 进行蒙哥马利变换
        for (size_t j = 0; j < S; j++) {
            C += (uint64_t)t[j] + m * (uint64_t)n_temp[j];
            t[j] = DATA(C);
            C >>= W;
        }
        C += t[S];
        t[S] = DATA(C);;
        C >>= W;
        t[S + 1] += C;

        // 左移 t 数组，准备下一�?
        for (size_t j = 0; j < S + 1; j++) {
            t[j] = t[j + 1];
        }
        //t++;
    }

    //printf("T:\n");
    sm2_uint256_t temp_R;
    uint32_t temp = t[S];
    for (size_t i = 0; i < S; i++) {
        temp_R[i] = t[S - i - 1];
        //printf("%08x ", T[i]);
    }
    if (temp || compare_value(temp_R, p) > 0) {
        sub(temp_R, Prime, temp_R);
    }
    convert_uint32_to_uint64(temp_R, R);
    //free(t_base);
}
void montgomery_multiplication_modn_64_32(uint64_t* A, uint64_t* B, uint64_t* p, uint64_t* R) {

    if (sm2_uint256t_is_zero(A) || sm2_uint256t_is_zero(B)) {
        sm2_uint256t_set_zero(R);
        return;
    }
    int S = 8;

    uint32_t t[10] = { 0 };
    uint32_t T[9] = { 0 };
    uint64_t C;
    //uint64_t C; // uint128
    //sm2_uint128_t C, res;
    uint64_t carry = 0;

    sm2_uint256_t A_temp, B_temp, n_temp, n_inv_temp, p_temp;
    sm2_uint256_64t neg_inv;


    bool flag = false;
    // 改为小端数据进行操作
    int cnt = 6;
    for (int i = 0; i < 4; i++) {
        // 小端存储，将64位数据拆分成2�?32位，并以小端格式存储
        A_temp[cnt] = (uint32_t)(A[i] & 0xFFFFFFFF);        // �? 32 �?
        A_temp[cnt + 1] = (uint32_t)((A[i] >> 32) & 0xFFFFFFFF); // �? 32 �?

        B_temp[cnt] = (uint32_t)(B[i] & 0xFFFFFFFF);        // �? 32 �?
        B_temp[cnt + 1] = (uint32_t)((B[i] >> 32) & 0xFFFFFFFF); // �? 32 �?

        n_temp[cnt] = (uint32_t)(p[i] & 0xFFFFFFFF);        // �? 32 �?
        n_temp[cnt + 1] = (uint32_t)((p[i] >> 32) & 0xFFFFFFFF); // �? 32 �?

        n_inv_temp[cnt] = (uint32_t)(N_NEG_INV_64[i] & 0xFFFFFFFF);    // �? 32 �?
        n_inv_temp[cnt + 1] = (uint32_t)((N_NEG_INV_64[i] >> 32) & 0xFFFFFFFF); // �? 32 �?

        cnt -= 2; // 由于每个64位数拆分为两�?32位数，所以减�?2
    }
    //print_uint32_array(A_temp);
    //printf("\n");
    // 外层循环
    for (size_t i = 0; i < S; i++) {
        C = 0;
        // 内层循环：计�? A[j] * B[i] 并累加到 t[j]
        for (size_t j = 0; j < S; j++) {
            C = (uint64_t)t[j] + (uint64_t)A_temp[j] * (uint64_t)B_temp[i] + C;
            t[j] = DATA(C);
            C >>= W;
        }

        C += (uint64_t)t[S];
        //t[S] = C & ((1 << W) - 1);
        //t[S + 1] = C >> W;
        t[S] = DATA(C);
        t[S + 1] = C >> W;
        C = 0;

        // 计算 m = t[0] * n_inv (mod 2^W)
        uint32_t m = (t[0] * n_inv_temp[0]) & 0xffffffff;

        // 进行蒙哥马利变换
        for (size_t j = 0; j < S; j++) {
            C += (uint64_t)t[j] + (uint64_t)m * (uint64_t)n_temp[j];
            t[j] = DATA(C);
            C >>= W;
        }
        C += t[S];
        t[S] = DATA(C);;
        C >>= W;
        t[S + 1] += C;

        // 左移 t 数组，准备下一�?
        for (size_t j = 0; j < S + 1; j++) {
            t[j] = t[j + 1];
        }
    }

    //printf("T:\n");
    sm2_uint256_t temp_R;
    uint32_t temp = t[S];
    for (size_t i = 0; i < S; i++) {
        temp_R[i] = t[S - i - 1];
        //printf("%08x ", T[i]);
    }
    if (temp || compare_value(temp_R, p) > 0) {
        //sub(temp_R, Prime, temp_R);
        sub(temp_R, param_N, temp_R);
    }
    convert_uint32_to_uint64(temp_R, R);
}


/**
  * Multiplication in Montgomery Field(P)
  *
  * @param a 256 bit number
  * @param b 256 bit number
  * @param p Prime numbers in Montgomery domain(P)
  * @param result multiplying two 256 bit numbers and taking mold yields 256 bit numbers
  */


void mod_sqr_modp_64(sm2_uint256_64t a, sm2_uint256_64t p, sm2_uint256_64t result)
{
    //mod_mul(a, a, p, result);
    montgomery_multiplication_modp_64_32(a, a, p, result);
}
void mod_sqr_modn_64(sm2_uint256_64t a, sm2_uint256_64t p, sm2_uint256_64t result)
{
    //mod_mul(a, a, p, result);
    montgomery_multiplication_modn_64_32(a, a, p, result);
}

void mod_fast_exp_modp_64(sm2_uint256_64t a, uint64_t b, sm2_uint256_64t p, sm2_uint256_64t result) {
    sm2_uint256_64t base, res;
    // 初始�?
    memcpy(base, a, sizeof(sm2_uint256_64t));  // base = a
    (void)ILib_memcpy(res, N_ONE_64, sizeof(N_ONE_64));
    //num2mgmy(N_256_1, p, res);

    while (b > 0) {
        if (b & 1) {
            //mod_mul(res, base, p, res);
            montgomery_multiplication_modp_64_32(res, base, p, res);
        }
        //mod_mul(base, base, p, base);  // base = base^2
        montgomery_multiplication_modp_64_32(base, base, p, base);
        b = b >> 1;                // exp = exp / 2
    }
    memcpy(result, res, sizeof(sm2_uint256_64t));
}
void mod_fast_exp_modn_64(sm2_uint256_64t a, uint64_t b, sm2_uint256_64t p, sm2_uint256_64t result) {
    sm2_uint256_64t base, res;
    // 初始�?
    //memcpy(base, a, sizeof(sm2_uint256_64t));  // base = a
    sm2_uint256t_copy(a, base);
    (void)ILib_memcpy(res, N_ONE_MOD_n_64, sizeof(N_ONE_64));

    //num2mgmy(N_256_1, p, res);

    while (b > 0) {
        if (b & 1) {
            //mod_mul(res, base, p, res);
            montgomery_multiplication_modn_64_32(res, base, p, res);
        }
        //mod_mul(base, base, p, base);  // base = base^2
        montgomery_multiplication_modn_64_32(base, base, p, base);
        b = b >> 1;                // exp = exp / 2
    }
    //memcpy(result, res, sizeof(sm2_uint256_64t));
    sm2_uint256t_copy(res, result);
}


/**
* Model inversion operation Montgomery Domain(P)
*/

void mod_inverse_modn_64(sm2_uint256_64t a, sm2_uint256_64t p, sm2_uint256_64t result)
{
    sm2_uint256_64t p_2;

    sm2_uint256_64t temp;
    sm2_uint256_64t a64_1;
    sm2_uint256_64t a64_2;
    sm2_uint256_64t a64_3;
    sm2_uint256_64t a64_4;
    sm2_uint256_64t a64_5;
    sm2_uint256_64t a64_6;
    sm2_uint256_64t a64_7;
    sm2_uint256_64t subnum2 = {
        0, 0, 0, 0x2
    };
    sub_64(p, subnum2, p_2);
    // 计算 a^(2^64)，原先是分多步计算到a^(2^32)后再逐步平方得到更高次幂，现在直接计算到a^(2^64)
    mod_sqr_modn_64(a, p, a64_1);         // a64_1  a^(2^64)
    for (int i = 0; i < 63; i++) {
        mod_sqr_modn_64(a64_1, p, a64_1);
    }

    // 分别先计算出a^(2^64*i)的结果（i = 2, 3, 4等）
    montgomery_multiplication_modn_64_32(a64_1, a64_1, p, a64_2);  // a64_2  a^(2^128)
    for (int i = 0; i < 63; i++) {
        mod_sqr_modn_64(a64_2, p, a64_2);
    }
    montgomery_multiplication_modn_64_32(a64_2, a64_2, p, a64_3);  // a64_3  a^(2^192)
    for (int i = 0; i < 63; i++) {
        mod_sqr_modn_64(a64_3, p, a64_3);
    }
    montgomery_multiplication_modn_64_32(a64_3, a64_3, p, a64_4);  // a64_4  a^(2^256)
    for (int i = 0; i < 63; i++) {
        mod_sqr_modn_64(a64_4, p, a64_4);
    }

    // 再把结果和a^(p - 2)[i]的结果做乘，这里的指数计算部分需要根据新的结构体定义调整
    mod_fast_exp_modn_64(a, p_2[3], p, temp);
    mod_fast_exp_modn_64(a64_1, p_2[2], p, result);
    montgomery_multiplication_modn_64_32(temp, result, p, temp);
    mod_fast_exp_modn_64(a64_2, p_2[1], p, result);
    montgomery_multiplication_modn_64_32(temp, result, p, temp);
    mod_fast_exp_modn_64(a64_3, p_2[0], p, result);
    montgomery_multiplication_modn_64_32(temp, result, p, result);
}
//void mod_inverse_64(sm2_uint256_64t a, sm2_uint256_64t p, sm2_uint256_64t result) {}
void mod_inverse_modp_64(sm2_uint256_64t in, sm2_uint256_64t p,sm2_uint256_64t out)
{
    sm2_uint256_64t t0, t1, t2, t3;
    sm2_uint256_64t ftmp;
    unsigned i;

    /* Step 1: t0 = a^3 = (2^2 - 2^0) * a */
    //felem_square(tmp, in);
    //felem_reduce(ftmp, tmp);
    //felem_mul(tmp, ftmp, in);
    //felem_reduce(t0, tmp);
    mod_sqr_modp_64(in, p, ftmp);
    montgomery_multiplication_modp_64_32(ftmp, in, p, t0);
    /* Step 2: t1 = t0^2 * a = (2^3 - 2^0) * a */
    //felem_square(tmp, t0);
    //felem_reduce(ftmp, tmp);
    //felem_mul(tmp, ftmp, in);
    //felem_reduce(t1, tmp);
    mod_sqr_modp_64(t0, p, ftmp);
    montgomery_multiplication_modp_64_32(ftmp, in, p, t1);
    /* Step 3: t2= t1^(2^3) * t1 = (2^6 - 2^0) * a */
    //felem_assign(ftmp, t1);
    (void)ILib_memcpy(ftmp, t1, sizeof(sm2_uint256_64t));
    for (i = 0; i < 3; i++) {
        //felem_square(tmp, ftmp);
        //felem_reduce(ftmp, tmp);
        mod_sqr_modp_64(ftmp, p, ftmp);
    }
    //felem_mul(tmp, ftmp, t1);
    //felem_reduce(t2, tmp);
    montgomery_multiplication_modp_64_32(ftmp, t1, p, t2);
    /* Step 4: t1 = t2^2 = (2^7 - 2^1) * a */
    //felem_square(tmp, t2);
    //felem_reduce(t1, tmp);
    mod_sqr_modp_64(t2, p, t1);
    /* Step 5: t3 = t1 * z = (2^7 - 2^0) * a */
    //felem_mul(tmp, t1, in);
    //felem_reduce(t3, tmp);
    montgomery_multiplication_modp_64_32(in, t1, p, t3);
    /* Step 6: t1= t1^(2^5) * t2 = (2^12 - 2^0) * a */
    //felem_assign(ftmp, t1);
    (void)ILib_memcpy(ftmp, t1, sizeof(sm2_uint256_64t));
    for (i = 0; i < 5; i++) {
        //felem_square(tmp, ftmp);
        //felem_reduce(ftmp, tmp);
        mod_sqr_modp_64(ftmp, p, ftmp);
    }
    //felem_mul(tmp, ftmp, t2);
    //felem_reduce(t1, tmp);
    montgomery_multiplication_modp_64_32(ftmp, t2, p, t1);
    /* Step 7: t2= t1^(2^12) * t1 = (2^24 - 2^0) * a */
    //felem_assign(ftmp, t1);
    (void)ILib_memcpy(ftmp, t1, sizeof(sm2_uint256_64t));
    for (i = 0; i < 12; i++) {
        //felem_square(tmp, ftmp);
        //felem_reduce(ftmp, tmp);
        mod_sqr_modp_64(ftmp, p, ftmp);
    }
    //felem_mul(tmp, ftmp, t1);
    //felem_reduce(t2, tmp);
    montgomery_multiplication_modp_64_32(ftmp, t1, p, t2);
    /* Step 8: t1= t2^(2^7) * t3 = (2^31 - 2^0) * a */
    //felem_assign(ftmp, t2);
    (void)ILib_memcpy(ftmp, t2, sizeof(sm2_uint256_64t));
    for (i = 0; i < 7; i++) {
        //felem_square(tmp, ftmp);
        //felem_reduce(ftmp, tmp);
        mod_sqr_modp_64(ftmp, p, ftmp);
    }
    //felem_mul(tmp, ftmp, t3);
    //felem_reduce(t1, tmp);
    montgomery_multiplication_modp_64_32(ftmp, t3, p, t1);
    /* Step 9: t2 = t1^4 = (2^33 - 2^2) * a */
    //felem_square(tmp, t1);
    //felem_reduce(ftmp, tmp);
    //felem_square(tmp, ftmp);
    //felem_reduce(t2, tmp);
    mod_sqr_modp_64(t1, p, ftmp);
    mod_sqr_modp_64(ftmp, p, t2);
    /* Step 10: t3= t2^(2^29) = (2^62 - 2^31) * a */
    //felem_assign(t3, t2);
    (void)ILib_memcpy(t3, t2, sizeof(sm2_uint256_64t));
    for (i = 0; i < 29; i++) {
        //felem_square(tmp, t3);
        //felem_reduce(t3, tmp);
        mod_sqr_modp_64(t3, p, t3);
    }
    /* Step 11: t1 = t1 * t3 = (2^62 - 2^0) * a */
    //felem_mul(tmp, t1, t3);
    //felem_reduce(t1, tmp);
    montgomery_multiplication_modp_64_32(t1, t3, p, t1);
    /* Step 12: t3 = t3^4 = (2^64 - 2^33) * a */
    //felem_square(tmp, t3);
    //felem_reduce(t3, tmp);
    //felem_square(tmp, t3);
    //felem_reduce(t3, tmp);
    mod_sqr_modp_64(t3, p, t3);
    mod_sqr_modp_64(t3, p, t3);
    /* Step 13: t0 = t0 * t3 * t2 = (2^64 - 2^0) * a */
    //felem_mul(tmp, t0, t3);
    //felem_reduce(t0, tmp);
    //felem_mul(tmp, t0, t2);
    //felem_reduce(t0, tmp);
    montgomery_multiplication_modp_64_32(t0, t3, p, t0);
    montgomery_multiplication_modp_64_32(t0, t2, p, t0);
    /* Step 14: t2= ((t3^(2^32) * t0)^(2^64) * t0)^(2^94) = (2^254 - 2^222 - 2^94) * a */
    //felem_assign(ftmp, t3);
    (void)ILib_memcpy(ftmp, t3, sizeof(sm2_uint256_64t));
    for (i = 0; i < 32; i++) {
        //felem_square(tmp, ftmp);
        //felem_reduce(ftmp, tmp);
        mod_sqr_modp_64(ftmp, p, ftmp);
    }
    //felem_mul(tmp, ftmp, t0);
    //felem_reduce(ftmp, tmp);
    montgomery_multiplication_modp_64_32(ftmp, t0, p, ftmp);
    for (i = 0; i < 64; i++) {
        //felem_square(tmp, ftmp);
        //felem_reduce(ftmp, tmp);
        mod_sqr_modp_64(ftmp, p, ftmp);
    }
    //felem_mul(tmp, ftmp, t0);
    //felem_reduce(ftmp, tmp);
    montgomery_multiplication_modp_64_32(ftmp, t0, p, ftmp);
    for (i = 0; i < 94; i++) {
        //felem_square(tmp, ftmp);
        //felem_reduce(ftmp, tmp);
        mod_sqr_modp_64(ftmp, p, ftmp);
    }
    //felem_assign(t2, ftmp);
    (void)ILib_memcpy(t2, ftmp, sizeof(sm2_uint256_64t));
    /* Step 15: out = (t1 * t2)^4 * a = (2^256 - 2^224 - 2^96 + 2^64 -1) * a */
    //felem_mul(tmp, t1, t2);
    //felem_reduce(ftmp, tmp);
    montgomery_multiplication_modp_64_32(t1, t2, p, ftmp);
    //felem_square(tmp, ftmp);
    //felem_reduce(ftmp, tmp);
    mod_sqr_modp_64(ftmp, p, ftmp);
    //felem_square(tmp, ftmp);
    //felem_reduce(ftmp, tmp);
    mod_sqr_modp_64(ftmp, p, ftmp);
    //felem_mul(tmp, ftmp, in);
    //felem_reduce(out, tmp);
    montgomery_multiplication_modp_64_32(ftmp, in, p, out);
}




/**
* Compare two uint32_t number
*/
bool compare(uint8_t* r, uint8_t* k, size_t len)
{
    if (r == NULL || k == NULL || len <= 0)
        return false;
    //return memcmp(r, k, len);
    return ILib_memcmp(r, k, len);
}

static uint64_t is_zero(uint64_t in)
{
    in |= (0 - in);
    in = ~in;
    in >>= 63;
    return in;
}

uint64_t sm2_uint256t_is_zero(const sm2_uint256_64t a)
{
    return
        is_zero(a[0]) &
        is_zero(a[1]) &
        is_zero(a[2]) &
        is_zero(a[3]);
}

void sm2_uint256t_copy(sm2_uint256_64t a, sm2_uint256_64t r)
{
    r[3] = a[3];
    r[2] = a[2];
    r[1] = a[1];
    r[0] = a[0];
}

void sm2_uint256t_set_zero(sm2_uint256_64t r)
{
    r[0] = 0;
    r[1] = 0;
    r[2] = 0;
    r[3] = 0;
}

bool is_in_N(uint64_t* r) {
    // 判断 r 是否 >= 1
    if (is_zero(r)) {
        return false;  // 如果 r 是零
    }

    // 判断 r 是否 < n
    if ((compare_value_64(param_N_64, r))) {
        return true;
    }

    return false;
}

bool is_point_on_curve(sm2_point_64t* point) {
    sm2_uint256_64t x, y, x_2, y_2, x_3, ax, tmp, y_2_tmp;

    //convert_uint64_to_uint8(point->x, x_bytes);
    //convert_uint64_to_uint8(point->y, y_bytes);
    num2mgmy_modp_64(point->x, Prime_64, x);
    num2mgmy_modp_64(point->y, Prime_64, y);

    // 2. 计算左侧：y^2
    mod_sqr_modp_64(y, Prime_64, y_2);

    // 3. 计算右侧：x^3 + ax + b
    mod_sqr_modp_64(x, Prime_64, x_2);
    montgomery_multiplication_modp_64_32(x_2, x, Prime_64, x_3);
    montgomery_multiplication_modp_64_32(x, mg_A_64, Prime_64, ax);
    mod_add_64(x_3, ax, Prime_64, tmp);
    mod_add_64(tmp, mg_B_64, Prime_64, y_2_tmp);

    // 4. 比较左侧和右侧是否相�?
    if (ILib_memcmp(y_2_tmp, y_2, sizeof(y_2)) == 0) {
        return true;  // 点在椭圆曲线�?
    }
    else {
        return false;  // 点不在椭圆曲线上
    }
}

void point_inverse_64(const sm2_point_64t* a, sm2_uint256_64t p, sm2_point_64t* r) {
    //mod_sub(p, a->x, p, r->x);
    (void)ILib_memcpy(r->x, a->x, sizeof(a->x));
    (void)ILib_memcpy(r->z, a->z, sizeof(a->z));
    mod_sub_64(p, a->y, p, r->y);
    //neg_p_64(r->y, a->y);
}

/**
 * Add two points
 *
 * @param a,b points in the Montgomery domain
 * @param p Prime numbers in Montgomery domain
 * @param result the result of add two points
 */


void Crypto_ECC_jacobi_to_affine_64(sm2_uint256_64t x, sm2_uint256_64t y, sm2_uint256_64t z, sm2_uint256_64t p, sm2_point_64t* r) {
    // 计算Z的模p逆元
    //sm2_uint256_t N_ONE = { 0,0,0,0,0,0,0,1 };
    sm2_uint256_64t Z2_inv, Z3_inv, z_squared, z_cubic;

    mod_sqr_modp_64(z, p, z_squared);   // Z1^2
    montgomery_multiplication_modp_64_32(z, z_squared, p, z_cubic);   // Z1^3

    mod_inverse_modp_64(z_squared, p, Z2_inv);
    mod_inverse_modp_64(z_cubic, p, Z3_inv);
    //_felem_inv(z_squared, p, Z2_inv);
    //_felem_inv(z_cubic, p, Z3_inv);


    // 计算仿射坐标x = X / Z (mod p)，y = Y / Z (mod p)
    montgomery_multiplication_modp_64_32(x, Z2_inv, p, r->x);  // x = X / Z^2
    montgomery_multiplication_modp_64_32(y, Z3_inv, p, r->y);  // y = Y / Z^3

}

void Crypto_ECC_jacc_point_double_64(const sm2_point_64t* a, sm2_uint256_64t p, sm2_point_64t* r) {
    sm2_uint256_64t lambda, x3, y3, z3, x1, y1, z1;
    sm2_uint256_64t tmp;
    if (sm2_uint256t_is_zero(a->x) && sm2_uint256t_is_zero(a->y) && sm2_uint256t_is_zero(a->z)) {
        sm2_uint256t_set_zero(r->x);
        sm2_uint256t_set_zero(r->y);
        sm2_uint256t_set_zero(r->z);

        return;
    }
    sm2_uint256t_copy(a->x, x1);
    sm2_uint256t_copy(a->y, y1);
    sm2_uint256t_copy(a->z, z1);
    sm2_uint256_64t x1_squared, z1_squared, y1_squared, x1_y1_squared, a_z1, z1_forth, y1_forth, _3_x1_squared, lambda1, lambda2, lambda3, lambda1_squared, _2_lambda2, lambda2_x3, first_part;
    mod_sqr_modp_64(y1, p, y1_squared);  // Y1^2
    mod_sqr_modp_64(x1, p, x1_squared);   // X1^2
    mod_sqr_modp_64(z1, p, z1_squared); // Z1^2

    mod_sqr_modp_64(z1_squared, p, z1_forth); // Z1^4
    //montgomery_multiplication_modp_64_32(N_THREE_64, x1_squared, p, _3_x1_squared); // 3 * X1^2 
    mod_add_64(x1_squared, x1_squared, p, tmp);
    mod_add_64(tmp, x1_squared, p, _3_x1_squared);
    montgomery_multiplication_modp_64_32(mg_A_64, z1_forth, p, a_z1); //a * Z1^4

    mod_add_64(_3_x1_squared, a_z1, p, lambda1);

    montgomery_multiplication_modp_64_32(x1, y1_squared, p, x1_y1_squared);
    
    mod_add_64(x1_y1_squared, x1_y1_squared, p, tmp);
    mod_add_64(tmp, tmp, p, lambda2);
    //montgomery_multiplication_modp_64_32(N_FOUR_64, x1_y1_squared, p, lambda2);// 4 * X1 * Y1^2

    mod_sqr_modp_64(y1_squared, p, y1_forth); // Y1^4
    //montgomery_multiplication_modp_64_32(N_EIGHT_64, y1_forth, p, lambda3); // 8 * Y1^4
    mod_add_64(y1_forth, y1_forth, p, tmp);
    mod_add_64(tmp, tmp, p, tmp);
    mod_add_64(tmp, tmp, p, lambda3);

    mod_sqr_modp_64(lambda1, p, lambda1_squared);
    //montgomery_multiplication_modp_64_32(N_TWO_64, lambda2, p, _2_lambda2); // 8 * Y1^4
    mod_add_64(lambda2, lambda2, p, _2_lambda2);
    mod_sub_64(lambda1_squared, _2_lambda2, p, x3);

    mod_sub_64(lambda2, x3, p, lambda2_x3);
    montgomery_multiplication_modp_64_32(lambda1, lambda2_x3, p, first_part);
    mod_sub_64(first_part, lambda3, p, y3);

    sm2_uint256_64t y1_z1;
    montgomery_multiplication_modp_64_32(y1, z1, p, y1_z1);
    mod_add_64(y1_z1, y1_z1, p, z3);

    //_jacobi_to_affine(x3, y3, z3, p, r);
    (void)ILib_memcpy(r->x, x3, sizeof(sm2_uint256_64t));
    (void)ILib_memcpy(r->y, y3, sizeof(sm2_uint256_64t));
    (void)ILib_memcpy(r->z, z3, sizeof(sm2_uint256_64t));

}

void Crypto_ECC_jacc_point_add_64(sm2_point_64t* a, sm2_point_64t* b, sm2_uint256_64t p, sm2_point_64t* r) {
    //sm2_uint256_t lambda, x1, y1, x2, y2, z1, z2 , x3, y3, z3;
    sm2_uint256_64t lambda, x3, y3, z3;
    if (sm2_uint256t_is_zero(a->x) || sm2_uint256t_is_zero(a->y)) {
        sm2_uint256t_copy(b->x, r->x);
        sm2_uint256t_copy(b->y, r->y);
        sm2_uint256t_copy(b->z, r->z);

        return;
    }
    if (sm2_uint256t_is_zero(b->x) || sm2_uint256t_is_zero(b->y)) {
        sm2_uint256t_copy(a->x, r->x);
        sm2_uint256t_copy(a->y, r->y);
        sm2_uint256t_copy(a->z, r->z);

        return;
    }

    //if (compare_value(x1, x2) == 0 && compare_value(y1, y2) == 0)
    if (ILib_memcmp(a->x, b->x, sizeof(sm2_uint256_64t)) == 0 && ILib_memcmp(a->y, b->y, sizeof(sm2_uint256_64t)) == 0)
    {
        //point_double(a, p, r);
        Crypto_ECC_jacc_point_double_64(a, p, r);
        return;
    }

    sm2_uint256_64t x1_squared, z1_squared, y1_squared, x1_y1_squared, a_z1, z1_forth, y1_forth, lambda1, lambda2, lambda3, lambda1_squared, _2_lambda2, lambda2_x3, first_part;
    sm2_uint256_64t z2_squared, z2_cubic, z1_cubic, lambda4, lambda5, lambda6, lambda7, lambda8, lambda9;
    sm2_uint256_64t lambda6_squared, lambda3_squared, lambda7_lambda3_2, _2_x3, lambda9_6, lambda3_cubic, lambda8_3, y3_, z3_, _2_inv;

    mod_sqr_modp_64(a->y, p, y1_squared);  // Y1^2
    mod_sqr_modp_64(a->x, p, x1_squared);   // X1^2
    mod_sqr_modp_64(a->z, p, z1_squared); // Z1^2
    mod_sqr_modp_64(b->z, p, z2_squared); // Z2^2
    montgomery_multiplication_modp_64_32(a->z, z1_squared, p, z1_cubic); // Z1^3
    montgomery_multiplication_modp_64_32(b->z, z2_squared, p, z2_cubic); // Z2^3

    montgomery_multiplication_modp_64_32(a->x, z2_squared, p, lambda1);
    montgomery_multiplication_modp_64_32(b->x, z1_squared, p, lambda2);
    mod_sub_64(lambda1, lambda2, p, lambda3);
    montgomery_multiplication_modp_64_32(a->y, z2_cubic, p, lambda4);
    montgomery_multiplication_modp_64_32(b->y, z1_cubic, p, lambda5);
    mod_sub_64(lambda4, lambda5, p, lambda6);

    mod_add_64(lambda1, lambda2, p, lambda7);
    mod_add_64(lambda4, lambda5, p, lambda8);

    mod_sqr_modp_64(lambda6, p, lambda6_squared);
    mod_sqr_modp_64(lambda3, p, lambda3_squared);
    montgomery_multiplication_modp_64_32(lambda7, lambda3_squared, p, lambda7_lambda3_2);
    mod_sub_64(lambda6_squared, lambda7_lambda3_2, p, x3);

    //montgomery_multiplication_modp_64_32(N_TWO_64, x3, p, _2_x3);
    mod_add_64(x3, x3, p, _2_x3);
    mod_sub_64(lambda7_lambda3_2, _2_x3, p, lambda9);
    montgomery_multiplication_modp_64_32(lambda9, lambda6, p, lambda9_6);
    mod_sqr_modp_64(lambda3, p, lambda3_squared);
    montgomery_multiplication_modp_64_32(lambda3_squared, lambda3, p, lambda3_cubic);
    montgomery_multiplication_modp_64_32(lambda8, lambda3_cubic, p, lambda8_3);
    mod_sub_64(lambda9_6, lambda8_3, p, y3_);

    montgomery_multiplication_modp_64_32(y3_, N_TWO_INV_64, p, y3);

    montgomery_multiplication_modp_64_32(a->z, b->z, p, z3_);
    montgomery_multiplication_modp_64_32(z3_, lambda3, p, z3);


    sm2_uint256t_copy(x3, r->x);
    sm2_uint256t_copy(y3, r->y);
    sm2_uint256t_copy(z3, r->z);
}



/**
 * Multiply a point by a number
 *
 * Note: K does not need to be converted to Montgomery domain
 *
 * @param P a point in the Montgomery domain
 * @param k integer in the natural domain
 * @param p Prime numbers in Montgomery domain
 * @param result the result of multiplying a point and a number
 */
void demo_point_multiply_64(sm2_point_64t* P, sm2_uint256_64t k, sm2_uint256_64t p, sm2_point_64t* R)
{
    // todo 可以优化
    sm2_point_64t temp_point;
    sm2_point_64t Q;
    //memset(&Q, 0, sizeof(sm2_point_t)); // 无穷远点 (0, 0)
    (void)ILib_memset(&Q, 0, sizeof(sm2_point_64t));
    (void)ILib_memset(&temp_point, 0, sizeof(sm2_point_64t));

    // 从最高位开始，逐位处理标量 k
    for (int i = 255; i >= 0; i--) {
        // 每次都做点的倍加
        Crypto_ECC_jacc_point_double_64(&Q, p, &Q);
        if (i == 255)
        {
            (void)ILib_memcpy(&Q, P, sizeof(sm2_point_64t));
        }
        // 判断标量当前位是否为1
        if (k[(i / 64)] >> (63 - (i % 64)) & 1) {
            // 如果当前位是1，就进行点加�?
            if (i == 255)
            {
                //Crypto_ECC_jacc_point_add_64(&ZERO_POINT_64, P, p, &temp_point);
                ILib_memcpy(temp_point.x, P->x, sizeof(sm2_uint256_64t));
                ILib_memcpy(temp_point.y, P->y, sizeof(sm2_uint256_64t));
                ILib_memcpy(temp_point.z, P->z, sizeof(sm2_uint256_64t));
            }
            else {
                Crypto_ECC_jacc_point_add_64(&Q, &temp_point, p, &temp_point);
            }

        }
    }

    // 返回结果
    memcpy(R, &temp_point, sizeof(sm2_point_64t));
}

void print_sm2_uint256_t(sm2_uint256_t num) {
    printf("random_k = { ");
    for (int i = 0; i < 8; i++) {
        printf("0x%08X", num[i]);
        if (i < 7) {
            printf(", ");
        }
    }
    printf(" }\n");
}

/**
* Print point coordinates
*/
void print_point(const sm2_point_t* Point)
{
    printf("    .x = ");
    printf("{");
    print_uint32_array(Point->x);
    printf("},\n");
    printf("\n    .y = ");
    printf("{");
    print_uint32_array(Point->y);
    printf("},\n");
    printf("\n    .z = ");
    printf("{");
    print_uint32_array(Point->z);
    printf("},\n");
    //printf("\n");
}

void print_point_64(sm2_point_64t* Point)
{
    printf("{    .x = ");
    printf("{");
    print_uint64_array(Point->x);
    printf("    },\n");
    printf("\n    .y = ");
    printf("{");
    print_uint64_array(Point->y);
    printf("    },\n");
    printf("\n    .z = ");
    printf("{");
    print_uint64_array(Point->z);
    printf("    }},\n");
    //printf("\n");
}

/**
* Print the hexadecimal representation of the uint32_t array which has 8 elements
*/
void print_uint32_array(uint32_t arr[8]) {
    for (int i = 0; i < 8; i++) {
        printf("0x%08x, ", arr[i]);
    }
    printf("\n");
}

void print_uint64_array(uint64_t arr[4]) {
    for (int i = 0; i < 4; i++) {
        printf( "0x%016llx, ", arr[i]);
    }
    printf("\n");
}

void print_uint128(uint64_t arr[2]) {
    for (int i = 0; i < 2; i++) {
        printf("%016llx ", arr[i]);
    }
    printf("\n");
}
/**
* Print the hexadecimal representation of the uint32_t array which has 16 elements
*/
void print_uint512_array(sm2_uint512_t arr) {
    for (int i = 0; i < 16; i++) {
        printf("%08x ", arr[i]);
    }
    printf("\n");
}

void print_uint512_array_64(sm2_uint512_64t arr) {
    for (int i = 0; i < 8; i++) {
        printf("%016llx ", arr[i]);
    }
    printf("\n");
}
/**
* Print the hexadecimal representation of the uint8_t array which has 32 elements
*/
void print_uint8_array(uint8_t arr[32]) {
    for (int i = 0; i < 32; i++) {
        printf("%02x ", arr[i]);
    }
    printf("\n");
}
#define CRYPTO_STOP_SEC_CODE
#include "Crypto_MemMap.h"
